Materials and method

2.1. Study area

The Western Ghats is one among the 35 global hotspots of biodiversity and it lies in the western part of peninsular India in a series of hills stretching over a distance of 1600 km from north to south and covering an area of about 1,60,000 sq.km. It harbors very rich flora and fauna and there are records of over 4600 species of flowering plants with 38% endemics, 330 butterflies with 11% endemics, 197 reptiles with 52% endemics, 529 birds with 5% endemics, 161 mammals with 12% endemics, 343 fishes with 39% endemics and 248 amphibians with 62% endemics (http://wgbis.ces.iisc.ernet.in/biodiversity/). Unplanned developmental activities during the post-industrialization period have threatened biodiversity, ecology and water availability. In this context, the government of India has created PAs network across wild habitats in order to conserve and reduce impacts. Fig. 2 depicts the locations of Kudremukh National Park (KNP), Rajiv Gandhi Tiger Reserve (RTR), Bandipur Tiger Reserve (BTR), which were chosen to understand the status of forests and LU dynamics. KNP comes under the Global Tiger Conservation Priority-I (as per Wildlife Conservation Society (WCS) and World Wide Fund-USA) and is a habitat to diverse endemic tax. RTR is also known as Nagarahole National Park, sensitive habitat to elephants (connects to Nilgiri Biosphere Reserve), and has been experiencing anthropogenic pressure, evident from degradation of forests, intensified agricultural activities. BTR is the first NP to get status as tiger reserve (TR) in 1973, which is also a habitat for numerous bird species (~250 species) and mammals endemic to Western Ghats. Rampant habitat degradations are common phenomenon across PAs due to unplanned developmental projects such as dams and interstate highways in the ecologically sensitive habitats apart from other pressures such as grazing, firewood collection and weed infestation (Eupatorium, Parthenium, Chromolaena odorata), illegal hunting, poaching (elephant poaching for ivory), etc. due to lack of effective ecomanagement leading to a very frequent instances of humananimal conflicts. Unauthorized occupation of forest lands for cultivation in the buffer region has aggravated human-animal conflicts leading to the frequent animal mortalities.

Fig. 2. Study area - select protected areas of Western Ghats, Karnataka.

2.2. Method

Mapping, monitoring, modelling and visualization of LULC dynamics in the select protected areas with 10 km buffer in the central Western Ghats has been carried out as discussed in Fig. 3. Remote sensing (RS) data used in the study are Landsat MSS (1973), TM (1992), OLI8 (2016), IRS p6L4X (2016), and online Google Earth data (http://earth.google.com). Topographic maps have provided ground control points to rectify RS data and scanned paper maps. Survey of India (SOI) toposheets (1:50000 and 1:250000 scales) and vegetation map of South India developed by French Institute (1986) of scale 1:250000 was digitized to identify various forest cover types and temporal analyses to find out the changes in vegetation. Pre-calibrated GPS (Global Positioning System - Garmin GPS unit) used for field measurements. Land use analyses involved
(i) generation of False Color Composite (FCC) of RS data (bandsegreen, red and NIR)
(ii) selection of training polygons by covering 15% of the study area (polygons are uniformly distributed over the entire study area)
(iii) loading these training polygons coordinates into pre-calibrated GPS,
(vi) collection of the corresponding attribute data (training) from field,
(iv) supplementing this information with Google Earth
(v) 60% of the training data has been used for supervised classification based on Gaussian maximum likelihood (GML) algorithm, while the balance is used for accuracy assessment.
GRASS GIS (Geographical Resources Analysis Support System, http://ces.iisc.ernet.in/grass), a free and open source software is used for analysis.

Fig. 3. Protocol for mapping, monitoring, modelling and visualization of LULC dynamics.

LST is estimated pixel-wise across varied land use categories as per the standard protocol (Bharath et al., 2013) based on each land use class emissivity that considers the amount of thermal infrared energy radiated from the Earth's surface of many different types of Earth surface processes, surface-atmosphere interactions. Estimation of LST using emissivity will provide more accurate estimation with appropriate calibration of atmospheric contamination by a separation of surface emissivity and temperature from radiance at ground level and atmospheric corrections (Li et al., 2013; Zak and O stir, 2012; Mohamed et al., 2017 Inamdar et al., 2008; Merlin et al., 2010). Emissivity is an intrinsic property of the surface, potential for enhancing ability to monitor landscape changes in environmentally sensitive zones as compared to other congruent practices (Hulley et al., 2014). The surface temperature of the area is calculated as per equation (3), through equations (1) and (2).

.............................(1)

.............................(2)

L_λ represents the top of atmosphere spectral radiance (Watts/ (m² srad μm)). λ represents the wavelength of the thermal band. Qcal is the DN value of the pixel under consideration. A_L and M_Lare offset and gain factor respectively whose value for TM is A_L = 1.238, M_L = 0.055 and OLI is A_L= 0.1, M_L = 3.3420 * 10^-4. K₁ (Watts/(m² srad mm)) and K₂ (kelvin) represents the pre-launch calibration coefficients. Value of K₁ is 607.76 and K₂ is 1260.56 (for TM), while K₁ for OLI (Band 11) is 480.8883 and K2 for is 1201.1442. TB represents the At-satellite brightness temperature. Emissivity correction is carried out using surface emissivity for the specified LC derived from the methodology described in (Stathopoulou et al., 2007. The emissivity corrected land surface temperature (Ts) is computed from equation (3).

.............................(3)

where, λ is the wavelength of emitted radiance for which the peak response and the average of the limiting wavelengths (l = 11. 5 mm) were used, ρ = h * c/σ (1.438 *10^-2 mK), σ = Stefan Boltzmann's constant (5.67 *10^-8 Wm^{- 2} K^-4 =1.38 * 10^-23J/K), h = Planck's constant (6.626 * 10^-34Jsec), c = velocity of light (2.998 * 10⁸m/sec), and ε is spectral emissivity.
CA is used (Equation (4)) to obtain a spatial context and distri- bution map is based on Markov transitional probability and area by combining multi- criteria land allocations to predict land use changes. The neighborhood effect is incorporated through diamond lter of a kernel size of 5 * 5 pixels to create spatially explicit contiguous weighing factors (Bharath et al., 2014; Chen et al., 2013; Grinblat et al., 2016; Matthews et al., 2007; Pontius and Malanson, 2005). The temporal land use analyses provided spatial and temporal aspects of LU changes and transition probability map and area matrix were generated based on Equations (5) and (6), Markov process. The original transition probability matrix (denoted by P) of land use type is obtained from two former land use maps (Grinblat et al., 2016; Matthews et al., 2007; Pontius and Malanson, 2005). CA coupled with Markov chain assisted in land use predictions of 2016 using the transitional probability and area matrix generated for 1973e1992. The validation has been carried out by evaluating with reference land use maps of 2016 (actual) using Kappa index for location and quantity. Based on these validations LU for 2025 was visualized considering intermediate iterations of 3-year time period. The CA model can be expressed as,

.............................(4)
where, S is the set of discrete cellular states, N is the Cellular field, t and t+1 indicate the different times, and F is the transformation function of cellular states in local space.
The Markov model is a theory based on the process of the formation of Markov random process systems for the prediction and optimal control theory method. Based on the Bayes conditional probability formula, the prediction of land use changes is calculated by the following equation:

.............................(5)

where, P_(N) is state probability of any times, and P _(N-1) is preliminary state probability. Transition area matrix can be obtained by,

.............................(6)

ij is the sum of areas from the i^th land use category to the j^th category during the years from start point to target simulation periods, and n is the number of land use types. The transition are a matrix must meet the following conditions